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Abstract. We discuss the physical basis of the statistical mechanics of self-gravitating systems. We show the cor- 
respondance between statistical mechanics methods based on the evaluation of the density of states and partition 
function and thermodynamical methods based on the maximization of a thermodynamical potential (entropy or 
free energy). We address the question of the thermodynamic limit of self-gravitating systems, the justification of 
the mean-field approximation, the validity of the saddle point approximation near the transition point, the lifetime 
of metastable states and the fiuctuations in isothermal spheres. In particular, we emphasize the tremendously long 
lifetime of metastable states of self-gravitating systems which increases exponentially with the number of parti- 
cles A'' except in the vicinity of the critical point. More specifically, using an adaptation of the Kramers formula 
justified by a kinetic theory, we show that the lifetime of a metastable state scales as e^^" in microcanonical 
ensemble and e'^^-' in canonical ensemble, where As and Aj are the barriers of entropy and free energy j = s — jdt 
(per particle) respectively. The physical caloric curve must take these metastable states (local entropy maxima) 
into account. As a result, it becomes multi- valued and leads to microcanonical phase transitions and "dinosaur's 
necks" (Chavanis 2002b, Chavanis & Rieutord 2003). The consideration of metastable states answers the critics 
raised by D.H.E. Gross |cond-mat/0 307535 /0403582]. 
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1. Introduction 

The statistical mechanics of self-gravitating systems has a 
. ^ long history starting with the seminal papers of Antonov 
^ , (1962) and Lynden-Bell & Wood (1968). A statistical me- 
$H ' chanics approach is particularly relevant to describe the 
. 5^ , late stages of "small" groups of stars {N ~ 10^), such as 
globular clusters, which evolve under the influence of stel- 
lar encounters ( "collisional" relaxation). Apart from as- 
trophysical applications, the statistical mechanics of stel- 
lar systems is of great interest in physics because it difi^ers 
in many respects from that of more familiar systems with 
short-range interactions (Padmanabhan 1990). In partic- 
ular, for systems with long-range interactions, the ther- 
modynamical ensembles are not equivalent, negative spe- 
cific heats are allowed in the microcanonical ensemble (but 
not in the canonical ensemble) and metastable equilibrium 
states can have tremendously long lifetimes making them 
of considerable interest. 

Two types of approaches have been developed to deter- 
mine the statistical equilibrium state of a self-gravitating 
system. In the thermodynamical approach, one determines 
the most probable distribution of particles by maximiz- 
ing the Boltzmann entropy at fixed mass and energy in 



the microcanonical ensemble or by minimizing the free 
energy F = E — TS at fixed mass and temperature in 
the canonical ensemble (Lynden-Bell & Wood 1968, Katz 
1978, Chavanis 2002a). This approach is the simplest and 
the most illuminating. In addition, it is directly related to 
kinetic theories (based on the Landau or on the Fokker- 
Planck equation) for which the Boltzmann entropy (or 
the Boltzmann free energy) plays the role of a Lyapunov 
functional and satisfies a H-theorem. Alternatively, in the 
statistical mechanics approach, one starts from the den- 
sity of states or partition function, transforms it into a 
functional integral and uses a saddle point approximation 
valid in a properly defined thermodynamic limit (Horwitz 
& Katz 1978, de Vega & Sanchez 2002, Katz 2003). 

In the first part of this paper, we discuss the connexion 
between these two procedures. We remain at a heuristic 
level, stressing more the physical ideas than the mathe- 
matical formalism. In Sec.|21 we introduce the entropy by 
a combinatorial analysis. In order to regularize the prob- 
lem at short distances, we consider either the case of self- 
gravitating fermions or the case of self-gravitating parti- 
cles with a soften potential. We also discuss the thermody- 
namic limit of the classical and quantum self-gravitating 
gas. In Sec. 13 we show the relation between the density of 
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states g{E) and the entropy functional S[f] and between 
the partition function Z(f3) and the free energy functional 
J[f] = S[f] — PE[f]. In the thermodynamic limit, the sad- 
dle point approximation amounts to maximizing the en- 
tropy at fixed mass and energy (microcanonical ensemble) 
or to minimizing the free energy at fixed mass and temper- 
ature (canonical ensemble). In Sec. ^ we discuss the no- 
tion of canonical and microcanonical phase transitions in 
self-gravitating systems. We perform the (standard) hori- 
zontal and (less standard) vertical Maxwell constructions 
and discuss the validity of the saddle point approxima- 
tion near the transition point for finite N systems. These 
results (e.g., microcanonical first order phase transitions) 
are relatively new in statistical mechanics and still subject 
to controversy (Gross 2003,2004). Therefore, we provide a 
relatively detailed discussion of these issues. 

In the second part of the paper, we emphasize the im- 
portance of metastable states in astrophysics and show 
how they can be taken into account in the statistical ap- 
proach. In Sec. 13 we use the Kramers formula to esti- 
mate the lifetime of a metastable state. We show that 
the lifetime of a metastable state scales as e^^'' in mi- 
crocanonical ensemble and e^^^ in canonical ensemble, 
where As and Aj are the barriers of entropy and free en- 
ergy j = s — (3e (per particle) respectively. Therefore, the 
typical lifetime of a metastable state scales as except 
in the vicinity of the critical point Ec (Antonov energy) 
or Tc (Emden- Jeans temperature). We explicitly compute 
the barriers of entropy and free energy close to the criti- 
cal point for classical self-gravitating particles (stars). The 
very long lifetime of metastable states, scaling as e^, was 
pointed out by Chavanis & Rieutord (2003) and the diffi- 
culty of a stellar system to overcome the entropic barrier 
and collapse was qualitatively discussed in Chavanis & 
Sommcria (1998). We here improve these arguments by 
developing a theory of fluctuations in isothermal spheres, 
following the approach of Katz & Okamoto (2000). We 
also determine how finite N effects affect the collapse tem- 
perature and the collapse energy. Finally, in Sec. we 
derive a Fokker-Planck equation for the evolution of the 
distribution of energies P(-B, t) in the canonical ensemble 
and make contact with the standard Kramers problem. 
We determine the typical lifetime of a metastable state by 
calculating the escape time accross a barrier of free energy. 

2. The most probable distribution 

2.1. The Fermi-Dirac distribution 

We consider a system of TV particles confined within a 
spherical box of radius R and interacting via Newtonian 
gravity. Let /(r,v, i) denote the distribution function of 
the system, i.e. f{r,v,t)cfir(Pv gives the mass of parti- 
cles whose position and velocity are in the cell (r, v; r -|- 
d^r, V + (Pv) at time t. The integral of / over the velocity 
determines the spatial density 



(1) 



and the total mass of the configuration is given by 

M = [ p d^r, (2) 



where the integral extends over the entire domain. On 
the other hand, in the meanfield approximation, the total 
energy of the system can be expressed as 



E 



fv'^d^rd'^v 



p^d^r = K + W, 



(3) 



where K is the kinetic energy and W the potential en- 
ergy. The meanfield expression of the potential energy is 
obtained from the exact expression 



W = (- 



1 



-GN{N- l)m^ 



1 



P2(ri,r2) 
|ri -rzl 



(4) 



by approximating the two-body distribution function 
-F2(ri,r2) by the product of two one-body distribution 
functions Pi(ri)xPi(r2) and using p{r) — NmPi. For 
self-gravitating systems, this mean-field approximation is 
exact in a proper thermodynamic limit N — > -|-oo with 
7? = l3GMm/R and A = -ER/GM'^ fixed (see Appendix 
IXji . The gravitational potential $ = —G J p(r')/|r— r'jd^r' 
is solution of the Newton-Poisson equation 



A<I> = AnGp. 



(5) 



In order to regularize the problem at short distances, 
we shall invoke quantum mechanics and use the Pauli 
exclusion principle. The Pauli exclusion principle is a 
fundamental concept in physics and it has also applica- 
tions in astrophysics, e.g. in white dwarf and neutron 
stars. Therefore, it can be considered as a physically rel- 
evant small-scale regularization for compact objects. We 
wish to determine the most probable distribution of self- 
gravitating fermions at statistical equilibrium. To that 
purpose, we divide the individual phase space {r, v} into a 
very large number of microcells with size [h/nif' where h 
is the Planck constant (the mass m of the particles arises 
because we use v instead of p as a phase space coordi- 
nate). A microcell is occupied either by or 1 fermion 
(or (7 = 2s + 1 fermions if we account for the spin). We 
shall now group these microcells into macrocells each of 
which contains many microcells but remains nevertheless 
small compared to the phase-space extension of the whole 
system. We call v the number of microcells in a macro- 
cell. Consider the configuration {ni} where there are ni 
fermions in the I*'' macrocell, n2 in the 2"*^ macrocell etc., 
each occupying one of the v microcells with no cohabita- 
tion. The number of ways of assigning a microcell to the 
first element of a macrocell is v, to the second v — 1 etc. 
Since the particles are indistinguishable, the number of 
ways of assigning microcells to all rii particles in a macro- 
cell is thus 

(6) 
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{v~ni)\ 
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To obtain the number of microstates corresponding to the 
macrostate {n-i} defined by the number of fermions rii in 
each macrocell (irrespective of their precise position in 
the cell), we need to take the product of terms such as 
^ over all macrocells. Thus, the number of microstates 
corresponding to the macrostate {ri^}, i.e. the probability 
of the state {rii}, is 



we find that the critical points of entropy correspond to 
the Fermi-Dirac distribution 



/ = 



1 + Ae^"('5-+*) ■ 



(15) 



w{{n,}) = n 



ni!(i/ - ni)! 



(7) 



This is the Fermi-Dirac statistics. As is customary, we 
define the entropy of the state {rii} by 



S{{n^}) = In W{{n,}). 



(8) 



It is convenient here to return to a representation in terms 
of the distribution function giving the phase-space density 
in the i-th macrocell 

-ii TO riirjo 



where A = e"''^ is a strictly positive constant (inverse fu- 
gacity) and /3 = -j^^ is the inverse temperature. Clearly, 
the distribution function satisfies f < 7]o, which is a con- 
sequence of Pauli's exclusion principle. 

So far, we have assumed that the system is isolated so 
that the energy is conserved. If now the system is in con- 
tact with a thermal bath (e.g., a radiation background) 
fixing the temperature, the statistical equilibrium state 
minimizes the free energy F = E — TS , or maximizes the 
Massieu function J = S — (3E, at fixed mass and temper- 
ature: 



Max J[f] I M[f]^M. 



(16) 



(A)3 



^g-j Introducing Lagrange multipliers and writing the varia- 
tional principle in the form 



where we have defined 770 = m^/h^, which represents the 
maximum value of / due to Pauli's exclusion principle. 
Now, using the Stirling formula, we have 



6 J 



0, 



(17) 



i^(ln u — 1) — vl — 



Vo 

ln<|j.( 1- A 







Kf) 


- 1 







Passing to the continuum limit ^ 0, we obtain the usual 
expression of the Fermi-Dirac entropy 



S^-kf 



1 - ^ 1 In 



we find that the critical points of free energy are again 
given by the Fermi-Dirac distribution H15() . Therefore, the 
critical points (first variations) of the variational prob- 
lems H1H(I and H16(l are the same. However, the stability of 
the system (regarding the second variations) can be differ- 

(10) 

ent in microcanonical and canonical ensembles (see, e.g., 
Chavanis 2002b). When this happens, we speak of a situa- 
tion of ensemble inequivalence. The stability of the system 
can be determined by a graphical construction, by simply 
(Prcfiv plotting the series of equilibria I3{E) and using the turn- 

ing point method of Katz (1978, 2003). Inequivalence of 

statistical ensembles occurs when the series of equilibria 
(11) presents turning points or bifurcations. 



h_\3 



(-) 



If we take into account the spin of the particles, the above 
expression remains valid but the maximum value of the 
distribution function is now 770 = gm'^/h^, where g — 2s+l 
is the spin multiplicity of the quantum states (the phase 
space element has also to be multiplied by g). In the non- 
degenerate (or classical) limit / <C r]Q, the Fermi-Dirac 
entropy Hll|l reduces to the Boltzmann entropy 



S = -kB I - 



In 



gm'^ 



1 



(12) 



2.2. Classical particles with soften gravitational 
potential 

We now consider a system of classical self-gravitating par- 
ticles, like stars in globular clusters. In order to make the 
problem of statistical mechanics well-posed mathemati- 
cally (see below), we introduce a soften potential of the 
form 



Now that the entropy has been precisely justified, the 
statistical equilibrium state (most probable state) of self- 
gravitating fermions is obtained by maximizing the Fermi- 
Dirac entropy (|lll) at fixed mass (|2J) and energy ©: 



y/{r - r')2 + rj 



(18) 



Max S[f] I E[f] - E, M[f] = M. 



(13) 



Introducing Lagrange multipliers 1/T (inverse temper- 
ature) and /i (chemical potential) to satisfy these con- 
straints, and writing the variational principle in the form 



6S-^6E+^SN^0, 



(14) 



where rg is the soften radius. As we shall see, the soften 
radius plays a role similar to the inverse of rjo, the max- 
imum phase space density, in the case of self-gravitating 
fermions. As said previously, this soften radius is intro- 
duced in order to pose the problem correctly. However, 
we shall argue in the sequel that this small-scale cut-off is 
irrelevant for the structure of stellar systems. 

We wish to determine the most probable distribution 
of stars at statistical equilibrium (Ogorodnikov 1965). To 
that purpose, we divide the individual phase space {r, v} 
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into a very large number of microcells with size {h/mf' 
where h is a constant with dimension of angular momen- 
tum. Of course, quantum mechanics is not relevant for 
stellar systems so that h should not be confused with the 
Planck constant in the present context. For classical sys- 
tems, a microcell can be occupied by an arbitrary number 
of particles. Adapting the counting analysis of Sec. 12.11 
to the present context, the number of microstates corre- 
sponding to the macrostate {ni}, i.e. its probability, is 



W{{ni\)^N\\{ 



nA 



(19) 



This is the Maxwell-Boltzmann statistics. If we define the 
entropy of the state {rii} by 



S{{n,}) = In W{{n, }), 



(20) 



and take the continuum limit, we obtain the usual expres- 
sion of the Boltzmann entropy 



/ 



In 



(21) 



Note that it differs from the expression H12|l obtained 
from the Fermi-Dirac entropy. This is of course related to 
the Gibbs paradox in standard thermodynamics (Huang 
1963). In the absence of self-gravity, Eq. H21(l reduces to 
the awkward expression 



S = NkBh-i 



V 



47rm E 



3/2- 



2 ' 



(22) 



which is clearly non-extensive. By constrast, Eq. (|12(l leads 
to the Sackur- Tetrode formula 



NkB In 



V fAnm E 



N\ 3h^ N 



3/2-1 



(23) 



which is extensive. As is well-known, the origin of this dis- 
crepency is due to the indiscernability of the particles and 
to the presence of the factor A^! in the Maxwell statistics 
(|19() . For a molecular gas, the Gibbs paradox is usualy 
solved by invoking quantum mechanics. For a system of 
stars, one cannot use this argument. We shall consider 
that the stars are discernable and use the expression (|21|l 
for the entropy. However, this choice does not affect the 
structure of the equilibrium state as we shall see in the 
sequel. 

The most probable distribution of stars at statistical 
equilibrium is now obtained by maximizing the Boltzmann 
entropy (|21|l at fixed mass and energy. This yields the 
Maxwell-Boltzmann distribution 



where $ is related to the density p by 



-G 



p(r') 



(24) 



(25) 



The microcanonical ensemble is the correct description of 
stellar systems which form an isolated Hamiltonian system 



in a first approximation. We can also consider the case of 
self-gravitating systems in contact with a thermal bath of 
non-gravitational origin which imposes its temperature T . 
For such systems, the correct description is the canonical 
ensemble and the statistical equilibrium state is obtained 
by minimizing the Boltzmann free energy F — E ~ TS 
at fixed mass. The canonical ensemble is also the correct 
description of a gas of self-gravitating Brownian particles 
(Chavanis, Rosier & Sire 2002). In this model, the friction 
and the stochastic fluctuations can mimick the influence 
of an external medium (thermostat) to which the system 
of origin is coupled. 

2.3. Thermodynamic limit of self-gravitating systems 

We introduce dimensionless variables such that r = Rr' , 
V = C/v' and / = {M/R^U^)f' where R is the box radius, 
M is the mass of the system and U = {GM/RY^^ is a 
typical velocity obtained by a Virial type argument (or 
dimensional analysis). For self-gravitating fermions, the 
entropy Hll|) can be expressed as 



S - 



-Nk 



dh'd^V 



x</:inf/: 

p. V M 



1 



/i 



In 1 



11. 



(26) 



where /i = ryoV G^MR^ is the degeneracy parameter 
(Chavanis & Sommeria 1998). Writing /i = {R/R^f^ 
with i?* — h? /GM^/^rrfi/^ ^ we note that the degeneracy 
parameter is the ratio, to the power 3/2, of the system's 
radius divided by the radius i?, of a "white dwarf star" 
(i.e. a completely degenerate ball of fermions) with mass 
M . The conservation of mass is equivalent to 



j f'd^r'dW 



= 1, 



(27) 



and the conservation of energy is equivalent to 



= / f — d'^r'd'^V 

GAf2 -'2 2 



P'{r'i)p'{r'2) 
|r'i-r'2| 



d^r[dh'^. 



Finally, the Massieu function can be written 

J^N{s[f]+rjA[f']), 



(28) 



(29) 



where s = S/N, t] = (3GMm/R and A = -ER/GAP. We 
define the thermodynamic limit as A^ — )■ +oo such that 
fj. = tjoVGHIW, a = -ER/GRP and 77 = PGMm/R are 
fixed. Coming back to physical quantities, it makes sense 
to fix h, m and G. Then, we have the scalings R ^ A^"^/"^, 
T - Af4/3, E - N^^ S N and J N a.s N ^ +00 
(the free energy F scales as N"^/"^). This is the quan- 
tum thermodynamic limit (QTL) for the self-gravitating 
gas (Chavanis 2002b, Chavanis & Rieutord 2003). This 
thermodynamic limit is relevant for compact objects with 
small radii R ^ N~^^^ <C 1 such as white dwarfs, neutron 
stars, fermion balls etc. The usual thermodynamic limit 
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N,R ^ +00 with N/ constant is clearly not relevant 
for inhomogeneous systems whose energy is non-additive 
(Padmanabhan 1990). 

For classical particles with soften potential, the en- 
tropy can be expressed as 



2.5 - 



S = -N f 



In 



- 1 



(30) 



where y = V G^MR^ is the counterpart of the degener- 
acy parameter. For classical particles, we see that it does 
not play any fundamental role in determining the struc- 
ture of the system since it just appears as an additional 
constant term (independent on /) in the entropy. If we 
only consider the part of entropy that depends on the dis- 
tribution function, we get 



Sr^-N j flnfd^r'd^V. 



(31) 



This is the relevant part of the entropy functional consid- 
ered by Antonov (1962) and Lynden-BeU & Wood (1968). 
We can therefore write S[f] = S'_r[/] + Sj where Sj is 
the constant part (irrelevant). The conservation of mass 
is equivalent to Eq. H27|l and the conservation of energy is 
equivalent to 



P'(r'i)p'(r'2) 



V'2? 



(32) 



where e — r^/R. As before, the Massieu function is 
given by Eq. (jJHJ- We define the thermodynamic limit 
as iV ^ +00 such that A = -ER/GM^, -q = (3GMm/R 
and e ~ tq/R are fixed. Coming back to physical quanti- 
ties, it makes sense to fix ro, m and G. Then, we have the 
scalings (Chavanis & Rieutord 2003) R ^ I, E ^ N^, 
T - A^, 5'i^ - TV and J^, - A^ as A^ -> +00 (the 
free energy Fr scales as A^^). These scalings imply that 
v ~ A^^/^ (if we fix h). Therefore, the (irrelevant) constant 
part of the entropy per particle diverges logarithmically as 
Si/N ^ \nv ^ i In A — ^ +00. This does not seem to be 
a crucial problem since this diverging term does not de- 
pend on / and therefore does not affects the structure of 
the equilibrium state. However, in a strict sense, there is 
no thermodynamic limit for classical self-gravitating par- 
ticles with soften potential. This contrasts with the case 
of self-gravitating fermions that possess a rigorous ther- 
modynamic limit (QTL). 

Let us finally consider the case of classical self- 
gravitating particles without small-scale cut-off. The en- 
tropy is given by the Boltzmann formula (|21|l . When vq — 
0, we know that the Boltzmann entropy has no global max- 
imum at fixed mass and energy (Antonov 1962). However, 
for sufficiently high energies, it has local entropy maxima 
that describe metastable gaseous states. The thermody- 
namic limit in that context corresponds to A — )■ +00 such 
that A = -ER/GAP and 77 = [iGMm/R are of order 
unity. If we fix TO, G and T, we have the scalings R ^ N, 



1.5 



0.5 











"gas 




Ay 


\ 


- 




LEM / 










/ SP 




- 


GEM^ 








Collapse 


Explosion 










E. 


LEM 






GEM 






~ ■ — 

C 




D 



-0.7 



-0.4 



-0.1 



0.2 



0.5 



Fig. 1. Series of equilibria for self-gravitating fermions 
with small cut-off/large /i/large system size R. It has a 
Z-shape structure (dinosaur's neck). There can be sev- 
eral values of inverse temperature (5 for a given energy E. 
They correspond to local maxima (LEM), global maxima 
(GEM) or saddle points (SP) of entropy S'[/]. The same 
remark applies in the canonical ensemble where the role 
of E and /? is reversed. 



5 - Af, 5 ~ A, J ~ TV and F - A as A ^ -Hoo. 
This is the classical thermodynamic limit (CTL), or di- 
lute limit, for the self-gravitating gas (de Vega & Sanchez 
2002). Physically, it describes metastable gaseous states 
that are not affected by the small-scale cut-off (Chavanis 

6 Rieutord 2003). As we shall see, these metastable states 
have considerably long lifetimes so that this thermody- 
namic limit is relevant for classical objects with large radii 
i? ^ A^ ^ 1 such as globular clusters. 

3. Connexion with statistical mechanics 

3.1. Series of equilibria and metastable states 

The critical points of entropy S[f] at fixed E and Af (i.e., 
the distribution functions /(r,v) which cancel the first 
order variations of S at fixed E, M) form a series of equi- 
libria parameterized, for example, by the density contrast 
TZ = p{0)/p{R) between the center and the edge of the 
system (see Chavanis 2002b). At each point in the series 
of equilibria corresponds a temperature (3 and an energy 
E. In this approach, (3 is the Lagrange multiplier asso- 
ciated with the conservation of energy in the variational 
problem (|14|l . It has also the interpretation of a kinetic 
temperature in the Fermi-Dirac distribution l|15|l . We can 
thus plot /3{E) along the series of equilibria. The form of 
this "caloric curve" depends on the value of the degener- 
acy parameter /i in the case of fermions (Chavanis 2002b) 
and on the soften radius e for regularized classical sys- 
tems (Chavanis & Ispolatov 2002). It also depends on the 
dimension of space D (Sire & Chavanis 2002, Chavanis 
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Fig. 2. Series of equilibria for self-gravitating fermions 
with large cut-off/small /i/small system size R. It has a 
A'^-shape structure. There is only one value of inverse tem- 
perature (3 for a given energy E. It corresponds to a global 
maximum of entropy (GEM). By contrast, there are sev- 
eral values of energy E for a given /3 in the canonical en- 
semble. They correspond to local maxima (LFEM), global 
maxima (GFEM) or saddle points (SP) of free energy J[/]. 
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Fig. 3. Series of equilibria for self-gravitating classical par- 
ticles without small-scale cut-off. There exists only lo- 
cal entropy maxima (metastable state) or unstable saddle 
points of entropy. The same remarks hold in the canonical 
ensemble. 



2004a). In D = 3, the caloric curve has a Z-shape (see 
Fig. ^1 for small cut-off and a A^-shape (see Fig. for 
large cut-off (for no cut-off, we recover the well-known 
spiral of Fig. There can be several values of temper- 
ature (3 for the same energy E because the variational 
problem H13|l can have several solutions: a local entropy 
maximum (metastable state) , a global entropy maximum. 



and one or several saddle points. We must represent all 
these solutions on the caloric curve because local entropy 
maxima (metastable states) are in general more physical 
than global entropy maxima for the timescales achieved 
in astrophysics. Indeed, the system can remain frozen in 
a metastable gaseous phase for a very long time. This is 
the case, in particular, for globular clusters and for the 
gaseous phase of fermionic matter (at high energy and 
high temperature) . The time required for a system placed 
in a metastable gaseous state to collapse is in general 
tremendously long and increases exponentially with the 
number N of particles. Thus, tji/e +oo in the thermo- 
dynamic limit N +00. The robustness of metastable 
states is due to the long-range nature of the gravita- 
tional potential. Therefore, at high temperatures and high 
energies, the global entropy maximum is not physically 
relevant. Condensed objects (e.g., planets, stars, white 
dwarfs, fermion balls,...) only form below a critical energy 
Ec (Antonov energy) or below a critical temperature Tc 
(Jeans temperature), when the gaseous metastable phase 
ceases to exist. The point where the metastable phase dis- 
appears is called a spinodal point. 

3.2. Microcanonical ensemble 

Let us explain things differently so as to make a closer 
contact with statistical mechanics. In statistical mechan- 
ics, one usually starts with the density of states 



N 

g{E)^ j J[i?-i/(ri,...,rw,vi,...,vjv)]n 



= 1 ^m' 



where H is the Hamiltonian. For our system 



H 



1 ^ 



(33) 



(34) 



4=1 



1<J 



The density of states is the normalization factor of the 
A^-body microcanonical distribution 



P/v(ri,Vi, ...,rAr,V7v) 



1 



9iE) 



SiE-H), 



(35) 



stating that all accessible microstates are equiprobable. 

Introducing the probability W{{ni}) of the state {n,}, 
we can rewrite the density of states in the form 



g{E) = Y W/({nJ), 

E{{n,})=E 



(36) 



where the sum runs over all macrostates with energy E. 
Introducing the entropy 5 = InW^ of the state {n^} and 
taking the continuum limit, the density of states can be 
expressed formally as 



g{E)^ / Vfe'^^f^SiE ~ E[f])S{M - M[f]), 



(37) 



where the sum runs over all distribution functions and 
S[f] is the Fermi-Dirac entropy Hll|l if the particles are 
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fermions and the Boltzmann entropy (|21|1 for classical 
particles (in that case, the gravitational potential must 
be regularized otherwise the density of states (|33|l di- 
verges). We now define the microcanonical entropy by 
S„iicro{E) ~ \ng{E) and the microcanonical temperature 
by Pmicro = dSmicro{E) / dE . By definition, the caloric 
curve (3micro{E) is uni- valued (Gross 2003). In the thermo- 
dynamic limit defined in Sec. 12.31 the entropy S[f] scales 
as ~ N, that is S[f] = Ns[f] where s 1 is the en- 
tropy per particle. Therefore, for N — s- -}-oo, the integral 
in Eq. (|S7|l is dominated by the state fgiobaii^i^) which 
is the global maximum of S[f] at fixed M and E (rig- 
orously speaking, we should work with the dimensionless 
quantities defined in Sec. 12.31 to get rid of the N — > +cx) 



dependence). Then, g(E) ~ e'^I-'^f'"''"'] , Sr, 



S[f global] 



and Pmicro = SS/SE = [3. However, this approach fails 
to take into account metastable states (local maxima of 
S[f] at fixed M and i?), which are of considerable interest 
in astrophysics. Indeed, equilibrium statistical mechanics 
tells nothing about timescales; a kinetic theory is required 
in that case. As explained above, these metastable states 
can persist for very long times. They correspond to the 
observed "diluted" structures in the universe (e.g., globu- 
lar clusters). Therefore, the caloric curve (3micro{E) does 
not describe the system adequately. The series of equi- 
libria 13{E) contain more information as they show local 
and global entropy maxima (as well as unstable saddle 
points). The curve f3micro{E) can be deduced from P{E) 
by keeping only global entropy maxima (see Fig. ^ . If we 
adopt this procedure, we find that the system exhibits a 
first order microcanonical phase transition (provided that 
the system size fi is sufficiently large) at a transition en- 
ergy Et{^) where the gaseous phase and the condensed 
phase have the same entropy (Chavanis 2002b). In the 
strict thermodynamic limit N — )■ +oo, this phase tran- 
sition is marked by a discontinuity of temperature. In 
fact, for finite N systems, the mean-field approximation 
g{E) ~ g^Slfgiobai] breaks down near the transition energy 
(Chavanis & Ispolatov 2002). This is because the contri- 
bution of the local entropy maximum in the functional 
integral H37(l becomes more and more important as we ap- 
proach Et- For the saddle point approximation to be valid, 
the number of particles must scale as N ~ |A — At|~^ 
for A ^ Af (see Sec. I4.2f) . Thus, for large but finite N, 
the temperature variation is sharp but remains continu- 
ous at the transition. We again emphasize that, due to 
the existence of metastable states, this phase transition 
may not be physically relevant. The true phase transition 
(gravothermal catastrophe) will rather take place at, or 
near, the spinodal point Ec (Antonov energy) where the 
metastable branch disappears. Estimating the lifetime of 
a metastable state by the Kramers formula tufe ^ e'^^' , 
where AS is the height of the entropic barrier (difference 
of entropy between the local maximum and the saddle 
point), we find that tuje ^ exp[2A'A^(Ac - A)^/^] with 
A' ~ 0.863159... (see SeclSJ. Except in the vicinity of the 
critical point Ac, the lifetime of a gaseous metastable state 
is enormous as it increases exponentially with the number 
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Fig. 4. Strict microcanonical caloric curve for ^ — 10^ 
(the control parameter is the energy E) . This figure is ob- 
tained from Fig. ^ by keeping only global entropy max- 
ima. It corresponds therefore to the true caloric curve 
(3micro{E) which is univalucd (Gross 2003). For N +00, 
there is a discontinuity of temperature at the transition 
energy Et{ii). For finite N systems, this discontinuity is 
smoothed out. Although this caloric curve is correct in a 
strict sense, it is not physical because it ignores metastable 
states that have an infinite lifetime in the thermodynamic 
limit. The physical caloric curve is obtained from Fig. ^ 
by discarding the unstable saddle points of entropy that 
form the intermediate branch (see Fig.lsj. 



of particles. Thus, metastable states have a considerable 
interest in astrophysics. 

If we now consider the case of classical particles {Ti ^ Q 
or /X — > +00), the transition energy Et{p) is rejected 
to +00 so that the whole branch of gaseous states is 
metastable. This "no cut-off" limit is relevant to classi- 
cal objects such as globular clusters or to the interstel- 
lar medium because, for these systems, the size of the 
particles (stars and atoms) clearly does not matter. In 
that case, the series of equilibria /3{E) forms a spiral (see 
Fig. |2Jl indicating the existence of one local entropy max- 
imum and one (or several) saddle points of entropy for 
a given energy (Lynden-Bell & Wood 1968). This spi- 
ral is the limiting form, for fi — > +00, of the fermionic 
caloric curve (see Fig. 11 in Chavanis 2002b). In this limit, 
the branch of "collapsed" states (condensed phase) co- 
incides with the X-axis where /3 = 0. It corresponds to 
configurations made of two particles in contact (^ bi- 
nary star) surrounded by a hot halo with T — > +00. 
This "binary+halo" configuration has an infinite entropy 
so, in a sense, it is the most probable configuration in 
the microcanonical ensemble (see Appendix A of Sire & 
Chavanis 2002). However, for sufficiently large energies 
(above the Antonov point), these configurations must be 
discarded "by hands" because they are reached for in- 
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Fig. 5. Physical microcanonical caloric curve for /i — 10^ 
(the control parameter is the energy E). This figure 
displays stable and metastable equilibrium states. The 
metastable states have infinite lifetime in the thermo- 
dynamic limit so they are physically relevant. The first 
order microcanonical phase transition of Fig. 0] does not 
take place in practice. Due to the existence of metastable 
states, the system displays a microcanonical hysteretic cy- 
cle marked by a "collapse" and an "explosion" at the spin- 
odal points where the branch of metastable states disa- 
pears (Chavanis & Rieutord 2003). 



accessibly large times. Therefore, for classical particles, 
the physical caloric curve Pphysicai (E) is obtained by tak- 
ing g{E) ~ e^^f'""''^ and Sphysicai ^ S[fiocai] where fiocai 
is the local entropy maximum at fixed mass and energy 
(Chavanis 2003). 
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Fig. 6. Strict canonical caloric curve for /i = 10'^ (the con- 
trol parameter is the inverse temperature (3). This figure 
is obtained from Fig. |2]by keeping only global maxima of 
free energy. It corresponds therefore to the true canonical 
caloric curve which is univalued and does not display neg- 
ative specific heats contrary to the corresponding micro- 
canonical caloric curve (see Fig. [Tjl. For N — > +oo, there 
is a discontinuity of energy (latent heat) at the transition 
temperature Tt{fj,). For finite N systems, this discontinu- 
ity is smoothed-out. Although this canonical caloric curve 
is correct in a strict sense, it is not physical because it 
ignores metastable states that have an infinite lifetime in 
the thermodynamic limit N — > +oo. The physical canon- 
ical caloric curve is obtained from Fig. [5] by discarding 
the unstable saddle points of free energy that form the 
intermediate branch with negative specific heats (see Fig. 

El). 



3.3. Canonical ensemble 

In the canonical ensemble, the object of fundamental im- 
portance is the partition function 



Z(^f3) ^ / g-/3H(ri,...,r„,vi,...,v„) 



DN 



(38) 



which is the normalization of the A^-body canonical dis- 
tribution 



PAr(ri, Vi, ...,rAr, Vat) 



1 



-I3H 



(39) 



The partition function can be rewritten 



(40) 



{".} 



where J{{ni}) = S{{ni}) — PE{{ni}) is the "free en- 
ergy" of the macrostate {rii} and the sum runs over all 



macrostates (in the present context, J[f] = S[f] — PE[f] 
is a more natural functional than the usual free energy 
F[f] ^ E[f] - TS[f]). Taking the continuum limit, the 
partition function can be expressed formally as 



ZiP)= J Vfe'^f^SiM ~ M[f]), 



(41) 



where the sum runs over all distribution functions and J[f] 
is the Fermi-Dirac free energy if the particles are fermions 
and the Boltzmann free energy for classical particles (in 
that case, the gravitational potential must be regularized 
otherwise the partition function (|38|l diverges). Note that 
Eq. H41|l can also be obtained from the exact formula 



Z{P) 



I e-0^giE)dE, 



(42) 



by substituting Eq. (|37|) for g{E) and carrying out the in- 
tegration over E. We now define the canonical free energy 
by Fcano — ~(l//3) i'D-Z. The average energy of the system 
at temperature T can be written {E)cano = —dhiZ/dp. 
By definition, the caloric curve {E)cano{P) is uni- valued. 
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Fig. 7. Microcanonical caloric curve for the same cut- 
off/degeneracy parameter ^/system's radius R as Fig. 
(the control parameter is the energy E) . It has a A^-shape 
structure and displays a stable region with negative spe- 
cific heats. All the equilibria are global maxima of en- 
tropy. In the canonical ensemble the region of negative 
specific heats is replaced by a phase transition (see Fig.O. 
The curves of Figs. and [3 are similar to those obtained 
by Padmanabhan (1990) with his toy model consisting in 
N = 2 stars in gravitational interaction. 
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Fig. 8. Physical canonical caloric curve for /i = 10'^ (the 
control parameter is the inverse temperature /3) . This fig- 
ure displays stable and metastable equilibrium states. The 
metastable states have infinite lifetime in the thermody- 
namic limit so they are physically relevant. The first order 
canonical phase transition of Fig. does not take place in 
practice. Due to the existence of metastable states the 
system displays a canonical hysteretic cycle marked by 
a "collapse" and an "explosion" at the spinodal points 
where the branch of metastable states disapears (Chavanis 
& Rieutord 2003). 



In the thermodynamic limit defined in Sec. 12.31 the free 
energy J[f] scales as ~ TV. Therefore, for N +oo, the 
integral in Eq. is dominated by the state /gZo&ai(r, v) 
which is the global maximum of J[f] at fixed M. Then, 

Z{P) ~ e-'I/.'"""-'], F,ano ^ -i^/P)J[f global] = E[f global] - 

TS[f global] and {E)cano = —SJ/SP = E. Metastable states 
(local maxima of J[f] at fixed M) can be taken into ac- 
count by plotting the full curve E{p). It is obtained from 
f3{E) defined in Sec. 13.11 by simply reversing the graph 
since the critical points of the variational problems 113|l 
and are the same (see Sec. l2.1|) . The curve {E)cano{P) 
can be deduced from E[I3) by keeping only global max- 
ima of free energy (see Fig. If we adopt this proce- 
dure, we find that the system exhibits a first order canon- 
ical phase transition at a transition temperature Tflp) 
where the gaseous phase and the condensed phase have 
the same free energy (Chavanis 2002b) . This phase transi- 
tion is marked by a discontinuity of energy (latent heat) . 
In fact, the mean-field approximation Z(P) ~ gJlfgiotai] 
breaks down near the transition temperature. For large 
but finite N, the energy variation is sharp but remains 
continuous at the transition. For the saddle point ap- 
proximation to be valid, the number of particles must 
scale as TV ~ \r] — r]t\^-^ for rj ^ r]t (see Sec. I4.2|) . We 
again emphasize that, due to the existence of metastable 
states, this phase transition may not be relevant and that 
the physical phase transition (isothermal collapse) takes 
place at, or near, the spinodal point Tc (Jeans temper- 
ature) where the metastable branch disappears (see Fig. 
IHl). Estimating the lifetime of a metastable state by the 
Kramers formula tufe g^F/ksT ^ where AF is the height 
of the potential barrier (difference of free energy between 
the local maximum and the saddle point), we find that 
tiije exp{2AA^(77c - rjf/^} with A ~ 0.16979815... (see 
Sec. EJ. Except in the vicinity of the critical point rjc, 
the lifetime of a gaseous metastable state is enormous as 
it increases exponentially with the number of particles. 
Metastable states are therefore highly robust. For classi- 
cal objects (?i 0), the transition temperature Tt(/i) is 
rejected to -l-oo so that the whole branch of gaseous states 
is metastable. In that case, the series of equilibria E{P) 
forms a spiral (see Fig. ^ indicating the existence of one 
local minimum of free energy F and one (or several) sad- 
dle points of free energy for a given temperature. In the 
classical limit, the branch of "collapsed" states (condensed 
phase) is rejected to E ^ —oo. It corresponds to config- 
urations where all the particles have collapsed at r = 0. 
This "Dirac peak" configuration has an infinite free energy 
F = —oo (due to the infinite binding energy) so, in a sense, 
it is the most probable configuration in the canonical en- 
semble (see Appendix B of Sire & Chavanis 2002). This 
differs from the binary star surrounded by a hot halo in the 
microcanonical ensemble. However, for sufficiently large 
temperatures (above the Emden- Jeans point), these con- 
figurations must be discarded "by hands" because they are 
reached for inaccessibly large times. Therefore, for classi- 
cal particles, the physical caloric curve EphysicaiiP) is ob- 
tained by taking Z{/3) ~ e'^[-^''"=°'l and Jphysicai ^ J[fiocai] 



10 



P.H. Chavanis: On the lifetime of metastable states in self-gravitating systems 



where fiocai is the local maximum of free energy J at fixed 
mass and temperature (Chavanis 2003). 

3.4. Grand canonical ensemble 

In the grand canonical ensemble, the partition fimction is 



(43) 



N=0 



Using Eq. ||1T|| . we get 



+ OC „ 

N=0 



= / T^fe 



Jlf]+0f,N[f] _ 



= / -Dfe 



Glf] 



(44) 



where G[f] = J[f] + (3^N[f] is the grand potential. Of 
course, the expression H41(l for is correct only for 
N 1. However, the contribution of small N terms in 
the grand partition function is expected to be weak so 
that Eq. provides a good approximation of the series. 
Now, the grand potential G[f] scales as ~ A^o = 'pU^- 
Therefore, in the thermodynamic limit R — s- +oo with 
fixed /3, G (gravitational constant) and m, the partition 
function Zqc is dominated by the distribution / which 
maximizes the grand potential G[f] at fixed (3 and ^. 
This problem has been considered for classical particles 
m D = 2> (Chavanis 2003). We shall reserve for a future 
work the study of self-gravitating fermions in the grand 
canonical ensemble. 



4. First order microcanonical and canonical phase 
transitions 

4.1. Maxwell constructions and critical points 

The deformation of the caloric curve when we vary the 
degeneracy parameter /z/ system size R is represented in 
Figs. El and ^| Similar curves are obtained for a hard 
sphere gas or a soften potential (Chavanis & Ispolatov 
2002) instead of fermions. In that case, l/fj, plays the role 
of the cut-off radius a or soften radius tq. 

For ^ < ^cTP — 82.5, the curve (3{E) defining the se- 
ries of equilibria is monotonic, so there is no phase transi- 
tion. For jjL > jiCTP, the curve E{l3) is multivalued so that 
a canonical first order phase transition is expected. The 
temperature of transition in the canonical ensemble can 
be obtained by a Maxwell construction as for the famil- 
iar Van der Waals gas. The equal area Maxwell condition 
^1 — ^2 (see Fig. ^ can be expressed as 

Ec 

(fi-l3t)dE^0, (45) 

Ea 

where Ea is the energy of the gaseous phase and Eb the 
energy of the condensed phase at the transition tempera- 
ture Tt . Since dS = pdE, one has 



Introducing the free energy J = S — PE, we verify that 
the Maxwell construction is equivalent to the equality of 
the free energy of the two phases at the transition: 



Ja = Jc- 



(47) 



If we keep only global maxima of free energy as in Fig. 
El the winding branch has to be replaced by a horizontal 
plateau. We see on Fig. that the extent of the plateau 
decreases as fi decreases. At the canonical critical point 
fJ-CTP, the plateau disappears and the curve presents an 
inflexion point. 

For fi > hmtp — 2600, the curve (}{E) is multival- 
ued so that a microcanonical first order phase transition 
is expected to occur (in addition to the canonical first or- 
der phase transition described previously). The energy of 
transition can be obtained by a vertical Maxwell construc- 
tion. The equal area Maxwell condition Ai = A2 (see Fig. 
[Tn|l can be expressed as 



/ {E - Et)dp = 0, 

•I Pa 



(48) 



where Ta and Tc are the temperatures of the two phases 
at the transition energy Et. Since dJ = —Edj3, one has 



Jc~JA + Et{l3c-f3A)^Q- 



(49) 



Thus, the Maxwell construction is equivalent to the equal- 
ity of the entropy of the two phases at the transition: 



Sa = Sc- 



(50) 



If we keep only global maxima of entropy as in Fig.^ the 
winding branch has to be replaced by a vertical plateau. 
We see that the extent of the plateau decreases as /i de- 
creases. At the microcanonical critical point ^j-mtp, the 
plateau disappears and the curve presents an inflexion 
point. 

Therefore, for fi > ^mtp, we expect a microcanonical 
and a canonical first order phase transition, for ^ctp < 

< fiMTP only a canonical first order phase transi- 
tion and for /i < ^ctp no phase transition at all. We 
emphasize, however, that due to the presence of long- 
lived metastable states, the first order phase transitions 
and the plateau are not relevant for the timescales of in- 
terest (see Sec. jS}. Only the zeroth order phase transi- 
tions (gravothermal catastrophe and isothermal collapse) 
marked by a discontinuity of entropy or free energy are 
physically relevant. 

4.2. Validity of the saddle point approximation near 
the transition point 

In this section, we discuss the validity of the saddle point 
approximation near the transition point. In the canonical 
ensemble, the partition fimction can be written 



Sc-SA-ptiEc-EA) = 0. 



(46) 



ZiP) 



g{E)dE, 



(51) 
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Fig. 9. Horizontal Maxwell plateau associated with a 
canonical first order phase transition. 
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Fig. 10. Vertical Maxwell plateau associated with a mi- 
crocanonical first order phase transition. 



where g{E) is the density of states with energy E. 
Introducing the entropy S{E) = lng{E), we can rewrite 
the partition function as 



J e'^'^^-P^'dE = J e'^'^UE = j c^'^^^Ue, 



(52) 



where J{E) = S{E) - (3E is the "free energy". As ex- 
plained in Sec. |2| the equilibrium states correspond to 
maxima of J. The condition J'{E) = 0, leading to S'{E) — 
(3, determines a series of equilibria E{f3). We shall con- 
sider the case where the series of equilibria E{(3) has the 
A^-shape structure of Fig. ^2 The transition temperature 
Tt{fj,), is determined by a Maxwell construction (see Sec 



I4.1|) . From Eq. H52I) . the distribution of energies at tem- 
perature T is given by 



P{E) 



Z{P) 



(53) 



It has a bimodal structure as shown in Fig. El The en- 
ergies where P{E) is maximum (denoted Ei and E2) cor- 
respond to the stable (GFEM) and metastable (LFEM) 
states on Fig. 1111 The energy where P{E) is minimum 
(denoted i?,) corresponds to the unstable states (SP) on 
Fig.m For T > Tt, gaseous configurations (high energies) 
are more probable than condensed configurations (low en- 
ergies). This is the opposite for T < Tt. Note that the 
notion of "more probable" is delicate since the system 
can remain blocked in a metastable ( "less probable" ) state 
for very long time, making that state the physically most 
likely state. 

For N +00, the partition function can be approxi- 
mated by 



Z{P) 



^Nj{E2) 



(54) 



where Ei and E2 are the energies at which J{E) is max- 
imum. We now wish to obtain the strict caloric curve 
{E)cano{P) defined in Sec. El When T is not too close 
from the transition temperature and N — > -foo, we 
need just keep the contribution of the global maximum 
of free energy as explained previously. To investigate the 
situation close to the transition temperature, we rewrite 
the partition function (|54(l as 



Z{P) 



1 



^N(](E2)-3(Et)) 



(55) 



Now, close to the transition point, we have j{Ei) — 
j{E2) + {rj — rjt) , where A is a constant of order unity 
depending on fi (see Fig. ^2 dashed fine). Therefore, 



1 



(56) 



where A77 — rj ~ rjf Thus, the saddle point approximation 
is valid for N\ri — r]t\ S> 1. This requires increasing large 
values of N as we approach the transition temperature Tt. 

In the microcanonical ensemble, the density of states 
can be written 



9{E) 



Vfe^^f^5{E ^ E[f])5{M - M[f]). 



(57) 



We shall consider the case where the series of equilibria 
f3{E) has the S'-shape structure of Fig.^l The transition 
energy Et{fi), is determined by a Maxwell construction 
(see Sec 14. If) . For A^ +00, the density of states can be 
approximated by 



5(£;) = e^^yil+e^^[^^l, 



(58) 



where /i and /2 are the distributions functions corre- 
sponding to the global and local maxima of entropy (we 
shall denote by /* the distribution function corresponding 
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Fig. 11. Caloric curve of the self-gravitating Fermi gas 
with /i = 1000. In the canonical ensemble, the temperature 
of transition is determined by a Maxwell construction. For 
J] < rjt {T > Tt), the gaseous states are global maxima of 
free energy J[f] and the condensed states local maxima 
(see the dashed curve). The situation is reversed for rj > rjt 
{T<Tt). 



to the saddle point of entropy). The corresponding tem- 
peratures (3i, P2 and /?* form the series of equilibria P{E) 
in Fig. El We now wish to obtain the strict caloric curve 
Pmicro{E) defined in Sec.|31 When E is not too close from 
the transition energy Et and N — > -l-oo, we need just keep 
the contribution of the global maximum of entropy. To in- 
vestigate the situation close to the transition energy, we 
rewrite the density of states H58|l as 



9{E) 



^ gNs[f\] 



1 



(59) 



Now, close to the transition point, we have s[fi] — s[f2] + 
A'^(A— At), where A' is a constant of order unity depending 
on fi (see Fig. Therefore, 



9iE) 



1 



where AA = A — At. Thus, the saddle point approximation 
is valid for iV|A — At| ^ 1. This requires increasing large 
values of N as we approach the transition point Et- 

5. The persistence of metastable states 

5.1. Typical lifetime of a metastable state 

In this section, we estimate the lifetime of a metastable 
state by using an adaptation of the Kramers formula 
(Risken 1989). We start first by the canonical ensem- 
ble. Close to the critical temperature Tc, see Fig. ^2 the 
free energy F[f] has one global minimum Fgiobai (stable 
GFEM), one local minimum Fiocai (metastable LFEM) 
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Fig. 12. Distribution of energies P{E) cx exp(5(i;) - /?£:) 
in the canonical ensemble for different values of temper- 
ature (here fx = 1000). The entropy S[E) has been cal- 
culated in the mean-field approximation, using Eq. 
The distribution of energies has a bimodal structure char- 
acteristic of a first order canonical phase transition. For 
V — Vt — 1.052, the two phases have the same probability. 
For r] = 1 < ?]t, the gaseous phase is the "most probable" 
and for 7y = 1.1 > r^t, the condensed phase is the "most 
probable" (in a strict sense). 

and one saddle point Fgaddie (unstable SP). We call Eiocai 
the energy of the metastable equilibrium state and Esaddie 
the energy of the unstable saddle point. For a system ini- 
tially prepared at Eiocai , the probability of the energy E 
is P{E) ~ e-^F{E)-Fio.ai)/kBT ^ js^Q^^ ^jjg energy fiuc- 

tuations drive the system past Eg^ddie^ it will collapse. 
Therefore, the lifetime of the metastable state can be es- 
timated by tufe ~ l/P{Esaddle) Or 
^AF/fesT 



tufe 



(61) 



(60) iuf 



where /S.F — \Fsaddie ~ ^ioco/l is the barrier of potential 
appropriate to our problem. Noting that 



(62) 



we conclude that, except in the vicinity of the critical point 
Tc, the lifetime of a metastable state scales as exp(7V). 
Therefore, metastable states (LFEM) are extremely ro- 
bust in astrophysics and cannot be neglected, even if there 
exists states with lower free energy (GFEM). To investi- 
gate the situation close to the critical point we shall 
calculate the barrier of potential Aj for the classical self- 
gravitating gas (?i ^ 0). To that purpose, we use the re- 
sults derived in a preceding paper (Chavanis 2002a). 

We recall that the series of equilibria is parameterized 
by a = (47rG/3po)^/^i?, where po is the central density. 
Introducing the Milne variables 



(63) 



P.H. Chavanis: On the lifetime of metastable states in self-gravitating systems 



13 



2.5 



DC 
5 

O 1.5 

CO. 

I 



0.5 









^=6000 


- 






- 


- 






- 




Tl. 


) 
















\ 

S(£) 





-0.7 -0.4 -0.1 0.2 
£=ER/GM' 



0.5 



0.8 



Fig. 13. Caloric curve of the self-gravitating Fermi gas 
with /i = 10^. In the microcanonical ensemble, the energy 
of transition is determined by a Maxwell construction. For 
e > et, the gaseous states are global maxima of entropy 
S[f] and the condensed states are local maxima (see the 
dashed curve). The situation is reversed for e < et. 



and noting uq — u{a) and vq = v(a) their value at the 
edge of the confining box, the thermodynamical parame- 
ters of the self-gravitating gas are given by 



A = — - uo 



5* — >S'o 



1 



(64) 
(65) 

Inuo - In(uo'yo) + t^o + 2mo — 3, (66) 



where 

^=ln,+ lln.-ln2-i. (67) 

The free energy per particle j = s + r^A is given by 

1 3 
3 = - 2 li^^o - In(uoi'o) + «o + uo - 2 ■ (^8) 

We now expand the Milne variables close to the crit- 
ical point Uc- Recalling that u(ct^ = 1, 'y(ac) = ''7c and 
v'{ac) — 0, we get 



uo = 1 + — (2 - ?7c)e - i^(2 - 

+ ^(2-r7c)(?7e+2r?c-4)e3 + .... 



(69) 



(70) 



where e — a — ac- We also recall that ac — 8.993195... and 
ric ~ 2.517551.... Substituting Eqs. ^ and ^ in Eq. 

we find that 



(71) 



where jo = ijc — (3/2)ln?7c — 1/2. On the other hand, 
recalling that rj ^ v{a), we get 



VciVc - 2) 2 VciVc - 4) 3 
Vc- V ^ e e 



2a? 



6ai? 



(72) 



Eliminating e from the foregoing relations, we finally ob- 
tain 



j = 30 + Vc) 



K ± A(7yc - v) 



1/2 



(73) 



with K = l/2?7c and A = (2/3)^2(?7c - '2)/r^l. We note 
that if = Aq where Aq is the normalized energy at 
the critical temperature (this comes from the fact that 
5J = —E5(3). Equation H73() reproduces the cusp at 77 = r/c 
formed by the curve j(ry) in Fig. 1111 Since dJ = —Edp, 
J(a) and 77(a) are extrema at the same points in the series 
of equilibria, which is at the origin of the cusp. From Eq. 
(|73|l . the barrier of free energy near the critical point r]c is 
given by 

Aj = 2A(r;,- 77)3/2. (74) 
Therefore, the lifetime of the metastable state scales as 



tlife 



(75) 



with A ~ 0.16979815.... Therefore, metastable states are 
robust if N{ric - ■qf'^ > 1. Note that if we had estimated 
the lifetime of the metastable state by tuf^ ~ g^E/ksT ^ 
we would have obtained tn/e ^ e"^^" with A" = 
[2(77c - 2)/77c]^/^ ^ 0.64121317.... We see that entropic 
effects modify the power oi rjc — rj in the expression of 
the metastable state lifetime. Returning to Eq. (TSJ, the 
temperature of collapse T; taking into account finite A'' 
effects can be estimated by 2XN{r]c — rji)^/^ ^ 1. This 
leads to 

2/3 



VI = Vc 



1 



1 

Vc 



N' 



-2/3 



A numerical application gives 
m = 2.517 (1 - 0.816 N-^/^). 



(76) 
(77) 



In the microcanonical ensemble, the lifetime of a 
metastable state can be estimated by 



tlife 



NAs 



(78) 



where AS* is the entropic barrier. By performing a study 
similar to the previous one close to the turning point of 
energy Ac (see Fig. lT^ . we find that 



So + (Ac - A) 



77o±A'(Ac-A)i/2 



(79) 
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with So = -0.192962. 
0.863159.... Therefore, 



As . 
and 

tlife 



2A'(Ae - A)3/2, 



= 2.03085... and A' = 



(80) 



(81) 



The energy of collapse Ei taking into account finite N 
effects can be estimated by 2A'A^(Ac - A)^/^ ^ 1. This 
leads to 



A; = A, 



1 



1 

a: 



1 

2A' 



2/3 



Ar-2/3 



a numerical application gives 

A/ = 0.335 (1 - 2.077 iV"^/^). 

This corresponds to a density contrast 

Ui = 708.6 (1 - 6.014 N-^/^). 



(82) 



(83) 



(84) 



These results are similar to those found by Katz & 
Okamoto (2000) by analyzing the microcanonical fluctua- 
tions of isothermal spheres. In particular, the scaling with 
N is the same. In the following section, we apply their 
theory of fluctuations to the canonical ensemble and show 
the relation with the preceding approach. 



5.2. Canonical fluctuations in isothermal spheres 

The canonical partition function can be written 

Z{/3) = / e'^'^^UE, 



(85) 



where J{E) = S{E) — (3E. As before, we shall consider 
the situation where J{E) has two maxima (stable) and one 
minimum (unstable). This corresponds to the caloric curve 
of Fig. ^] We shall be interested here by the metastable 
state, i.e. the /ocaZ maximum of J{E). Thus, we eliminate 
"by hands" the contribution of the global maximum in the 
integral. Expanding J{E) around the local maximum up 
to second order, we obtain 



,^j"(E)(SEf 



d{SE), 



(86) 



where E now refers to the energy of the metastable state 
and 5E is a fluctuation around equilibrium. According to 
the foregoing equation, the probability of the fluctuation 
SE is given by 



Noting that J"{E) = S"{E) and recalling that || = (3{E) 
at equilibrium, we can rewrite the foregoing expression as 



P{dE) = 



1 



/2^ 



dp 



dE 



1/2 



This formula shows that only equilibrium states with pos- 
itive specific heats C = dE/dT = —(3'^^ > are stable in 



the canonical ensemble. The variance of the fluctuations 
of energy is given by the usual formula 



1 



{{SEY) = - =t'c>o 



dE 



(89) 



which can also be directly derived from the exact A'^- 
body distribution function 139|) . These results are valid 
on the whole gaseous branch of Fig. ^2 If we now exam- 
ine the situation close to the critical point (/3c, Eq) where 
{d(3/dE)c = 0, we have to first order 



dp 




d^P 


dE 




dE^ 



[E - Eo) 



(90) 



We note E' the energy of the unstable saddle point of 
free energy J at temperature T. Close to the critical point 
{Pc Eq) we can approximate the caloric curve f3{E) by a 
parabole. Thus, E' is related to E, the energy of the local 
maximum of free energy J at temperature T, by 



E-E' = 2{E - Eo). 



(91) 



Using the criterion of Katz & Okamoto (2000), adapted 
to the canonical ensemble, we define the temperature of 
collapse as the temperature Ti at which the typical fluctu- 
ations of energy SE ~ {(SE)'^) are of the same order as 
the difference E — E' . Indeed, as we approach the critical 
point (PciEq) the fluctuations of energy become so im- 
portant (since the specific heat diverges) that the system 
can overcome the barrier of potential played by the saddle 
point of free energy and collapse (eventually reaching the 
global maximum of J). Thus, for finite N systems, grav- 
itational collapse can occur before the ending {/3c, Eq) of 
the metastable branch (spinodal point). According to the 
preceding criterion, the temperature of collapse is deter- 
mined by the condition 



{{6E)^)p,=4{Ei-Eo)'. 

Using Eqs. JBHl and (jHOJ, we obtain 



Ei=Eo 



d^(3 



dE^ 



-1/3 



(92) 



(93) 



It is more logical to write this criterion in terms of the 
temperature Ti. Close to the critical point, we have 



(3c -13 = 



{E - Eo)' 



dE^ 

Thus, using Eq. we get 

1/3 



(87) pi=p^-tU 



dE^ 



(94) 



(95) 



We can obtain more explicit results in the case of 
classical isothermal spheres (see Fig. ll4|) . Introducing the 
usual dimensionless parameters A = —ER/GM^ and 
T] = (3GMm/ R, the foregoing relation can be rewritten 

1/3 

(96) 



1 



8iV2/3 



d^TJ 



dA2 
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Fig. 14. Caloric curve showing metastable states (local 
maxima of J) and unstable equilibria (saddle points of J) 
in the canonical ensemble for classical particles {Ti ^ 0). 
As we approach the critical temperature r]c, the fluctua- 
tions of energy increase considerably allowing the system 
to collapse before the end of the metastable branch (spin- 
odal point). The temperature of collapse T; is the temper- 
ature at which the typical fluctuations of energy are of the 
same order as the energy difference |e' — e|. 



Now, we have (see Chavanis 2002a) 



(97) 



Therefore, the temperature of collapse for finite N systems 
is given by 

1/3-1 



1 



1 



(98) 



87V2/3 Vr?c-2^ 
A numerical application gives 

r/i = 2.517 (1 - 0.247 A^~2/3-)^ ^gg) 

This estimate can be compared with Eq. (|77|l which has 
the same scaling with N. Of course, we should not give too 
much credit on the numerical factor in front of N~'^/^ since 
the criterion for determining Ti is essentially phenomeno- 
logical. The corresponding energy A; can be deduced from 
Eq. (|94|l and is given by 

1/3 



A/ = Ao 



1 



■He 



(100) 



It may also be of interest to determine the corresponding 
density contrast TZi. To that purpose, we start from the 
formula (see Chavanis 2002a) 

„2 



a 



wofo 



(101) 



and we use the expansion ()69|l and ()70f) of the Milne vari- 
ables. This yields 



1 H e; 

etc 



(102) 



where TZc — ct'^/Vc- Now, e/ is determined by Eqs. H72|l 
and (IHHll yielding 



4l/37Vl/3(^^ _ 2)2/3- 

Substituting this result in Eq. H102|) we finally obtain 



(103) 



11c 



[2(r/c - 2)]2/3 
or, with numerical values. 
Til = 32.1 (1 - 2.45 N-^l^) 



(104) 



(105) 

The above theory thus predicts that, for finite N sys- 
tems, the collapse should take place slightly before the 
canonical spinodal point rjc due to the enhancement of 
energy fluctuations as we approach this critical point. 
The Monte Carlo simulations of de Vega & Sanchez 
(2002) in the canonical ensemble (with N — 2000 par- 
ticles) reveal that the collapse indeed takes place before 
the critical point. However, the collapse occurs appar- 
ently at the point where the isothermal compressibility 
kt ~ —{llV){dV/dp)T diverges. This corresponds to an 
inverse normalized temperature rjT = 2.43450.... It is sen- 
sibly smaller than the value obtained from Eq. (|99|l with 
N = 2000. The same discrepency with the prediction of 
Katz & Okamoto (2000) is found in the microcanonical en- 
semble by de Vega & Sanchez (2002). These results seem 
to indicate that the higher collapse temperature and en- 
ergy found in Monte Carlo simulations are not due to finite 
N effects. They seem to be independent on N and corre- 
spond to other critical points that do not coincide with the 
spinodal point. We intend to perform independent Monte 
Carlo simulations to check these results. 



5.3. General expression of the potential barrier 

We can use the preceding approach to obtain a simple 
approximate expression of the potential barrier close to 
the critical point rjc- Consider a system at fixed tempera- 
ture T and denote by E its equilibrium energy such that 
J'{E) — (we consider here that E is the energy of the 
metastable state). For a fluctuation E ~\- SE, we have seen 
that the variation of free energy SJ = ^J"{E){6E)^ can 
be expressed as 

SJ^^l3"{Eo)iE-E'){SEf 



(106) 



where we have assumed that E is close to the critical point 
Eq and, as before, E' is the energy of the minimum of J 
at temperature T. We can use this expression to estimate 
the potential barrier AJ = Jiocai — Jsaddie- Thus, setting 
SE = E- E', we get 



AJ=^\P"iEo)\iE~E'f. 



(107) 



Using Eq. H94|l to express this relation in terms of the 
temperature, we finally obtain 



AJ = 



32 



\f3"{Eo)\ 



(/3c - P) 



3/2 



(approx.). 



(108) 
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This is an estimate, because the curve J{E) is not just a 
parabole between E and E'. Using Eq. (|^ . this approx- 
imate expression can be written 



Aj = 6X{r]c - r]) 



3/2 



(approx.). 



(109) 



It differs from the exact expression H74|l by a factor 3. 
Noting that A J = ^{AE)^ / {{SEf) , according to Eqs. 

and (IHnil, the criterion of Katz & Okamoto 

(2000) corresponds to A J ~ 1/2. Alternatively, writing 
tiife ~ e^"^, the criterion leading to Eq. l(75)l corresponds 
to A J ~ 1. On a qualitative point of view, the approaches 
of Sees. 15. ll and 15.21 are equivalent. They slightly differ on 
the details (definition of collapse temperature, estimate of 
AJ) explaining why Eqs. (|77|l and H99|l are not exactly 
the same. 

We can also try to calculate AJ by working directly 
on the series of equilibria /3{E). Taking E a.s a con- 
trol parameter, we have J{E) = S{E) — [3{E)E. Noting 
that J'{E) = -(3'{E)E, J"{E) = -P"{E)E - /3'{E) and 
J"'{E) = -P"'{E)E-2f3"{E), and expanding J{E) close 
to the critical point where /3'{Eq) = 0, we get 



J{E) = JiEo) - -P"iEo)EoiE - Eo^ 
^^[f3"'{Eo)Eo + 2p"{Eo)]{E ~ Eo f + ... 



(110) 



Using the relation 

p-f3c = \p"iEo){E - E,f + ^^P"'{E^){E ~ E^f + ... 

(Ill) 

between the temperature and the energy close to the crit- 
ical point, we find that 



J{E) = J{Eo)-Eo{(3-P,) 



2(/3 - /3c 



3/2 



(112) 



Therefore, close to the critical point, the barrier of free 
energy is exactly given by 



AJ 



1 



32 



l/3"(i?o)| 



{(3c - Pf'^ (exact). 



(113) 



which returns Eq. (|74(l . 

We can use the same type of approach in the micro- 
canonical ensemble to obtain a simple approximate ex- 
pression of the entropic barrier close to the critical point 
(Ac, 770). Consider a system at fixed energy E and denote 
by (3 its equilibrium temperature (we consider here that 
(3 is the temperature of the metastable state). For a fluc- 
tuation (3 + 6(3, the variation of entropy can be expressed 
as (see Katz & Okamoto 2000 for details) 



5S^^E'{(3){S(3)\ 



(114) 



Assuming that /3 is close to the critical point /3o, and using 
arguments similar to those developed previously, we get 



6S^\E"iPo)iP-P')m' 



(115) 



where /3' is the inverse temperature of the saddle point 
of entropy. We can use this expression to estimate the 
entropic barrier AS* — Siocai ~ S saddle- Thus, setting S() — 
(3 — P', we get 



AS=\E"{(3o){(3~Pr- 
Using 



E~E,~-E"iPo){P-(3o)^ 



(116) 



(117) 



to express the temperature as a function of the energy 
close to the critical point, we finally obtain 



AS*: 



32 



{E - Ecf/"^ (approx.). 



(118) 



We note that Eqs. (|1()8|1 and (|118|l are symmetrical pro- 
vided that we interchange E and (3. Evaluating numeri- 
cally (PA/drf at the critical point {Ac,7]o), this approxi- 
mate expression can be written 



As = 6A'(A - A^f/^ (approx.). 



(119) 



It differs from the exact expression ()80|l by a factor 3. 

We can also try to calculate AS* directly from the series 
of equilibria E[(3). Taking /3 as a control parameter, we 
have S{P) = J(/3) + (3E{(3) and we recall that dS = (3dE 
and dJ = ^Ed(3. Thus S'{(3) = (3E'{(3), S"{(3) = E'{(3) + 
(3E"{(3) and S"'{(3) = 2E"{(3)+(3E"'{P). Expanding S{(3) 
close to the critical point where E'{(3q) = 0, we get 

5(/3) = 5(/3o) + ^£;"(/3o)/3o(/3-/3o)' 
+ |[2i?"(/3o) + /3o£;"'(/3o)](/3 - Pof + ... (120) 
Using the relation 

E-E, = \e"{Po){P - (3o)' + ^E"'i(3o)i(3 - (3o)' + ... 

(121) 

between the energy and the temperature close to the crit- 
ical point, we find that 



S{(3) = S{(3o) + (3o{E ~ E,) ± ^E"{(3o) 



2{E^E,) 



-1 3/2 



E"{(3o) 



(122) 



Therefore, close to the critical point, the barrier of entropy 
is exactly given by 



AS=l-, 



32 



Sy E"{(3o 
which returns Eq. (|80|l 



{E~Ecf/^ (exact), 



(123) 
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6. Relation to the Kramers problem 

6.1. The Fokker-Planck equation 

In the preceding section, we have used the Kramers for- 
mula to estimate the lifetime of metastable states in self- 
gravitating systems. We would like now to justify this for- 
mula from first principles. In order to determine the life- 
time of a metastable state, we need to introduce a dynam- 
ical model. In the canonical ensemble, we can consider a 
system of self-gravitating Brownian particles (Chavanis, 
Rosier & Sire 2002) described by the stochastic equations 

dvi 

(1^4) 

dvi _ 
dt 

where 



-V,C/(ri, ...,rA.) - Cv, + V2DR,{t), 



(125) 



-^Vi is a friction force and Ri(i) is a white noise 
satisfying (Ri(t)) = and {Ra,i{t)Rb,j{t')) = SijSabS{t ~ 
t'), where a,b = 1,2,3 refer to the coordinates of space and 
i,j = 1,...,N to the particles. The particles interact via 
the gravitational potential U (ri, tat) = J2i<j ^('"^i ~ ^j) 
where u{ri — Vj) = —G/\ri — rj\. The inverse temperature 
(3 = 1/T is related to the diffusion coefficient through the 
Einstein relation ^ = Dp. 

Using standard stochastic processes, we can derive the 
A^-body Fokker-Planck equation (Chavanis 2004b) 



9P, 



N 



dt 



N 

E 



dP, 



N 



N 



N 

E 



d 
dvi 



D 



.dPi 



N 



where F.^ 



dvi 

V,[/(r] 



(126) 



.,rAr) is the gravitational force 
acting on particle i and PAr(ri, vi, ...,rAr,VAr,t) is the 
-/V-body distribution function. Its stationary states cor- 
respond to the canonical distribution 



P7v(ri, vi, ...,rAr, Vat) 



1 



ZiP) 



(127) 

If we implement a mean-field approximation (Chavanis 
2004b), we can show that the distribution function 
f{r,v,t) = NPi is solution of the Kramers-Poisson sys- 
tem. However, this is not the approach that we shall con- 
sider here. 

We wish to obtain the time evolution of the distribu- 
tion of energies P{E,t). To that purpose, we shall follow 
the method developed by Kramers (1940) in his investiga- 
tion of the escape of Brownian particles over a potential 
barrier. The difference is that we work here in a 6N dimen- 
sional phase space. Assuming that P/v(ri, vi, rjv, vjv, <) 

depends only on energy E — J2iLi ^ + U{ri, r^r) and 
time t, and averaging the Fokker-Planck equation l|12t)|l 
on the hypersurface of energy i?, we show in Appendix B 
that 

OPn 



9{Ef-^{E,t) 



d 



m 



dE 



N 



where g{E) is the density of states and I{E) is the phase 
space hypervolume with energy less than E (thus g{E) = 
dl /dE). Now, the distribution of energies is given by 



P[E,t)=PN[E,t)g[E). 
At equilibrium, using Eq. (|127|l . we have 
1 



P{E) 



Z{(3) 



9{E)e 



-I3E 



(129) 



(130) 



Out of equilibrium, substituting Eq. I|129() into Eq. (|128f) 
and simplifying the resulting expressions, we finally obtain 



dP 
'dt 



_d_ 

'dE 



D{E)[^+(3PF'{E) 



(131) 



where D{E) = 3MI{E)/g{E) and F{E) = E- TS{E) = 
E — T In g{E) is the free energy. This is similar to the 
Fokker-Planck equation describing the stochastic motion 
of a particle in a potential where the energy E plays the 
role of the position x and where the free energy F{E) plays 
the role of the potential U{x). In the following, we shall 
assume that the free energy F{E) has a local minimum at 
Ea (metastable), a local maximum at Eb (unstable) and a 
global minimum at Eq. A typical situation is illustrated in 
Fig. El We shall prepare a large number A/" ^ 1 of systems 
close to the energy Ea with the canonical distribution 
(|130|l . Thus, AfxP{E,t)dE gives the number of systems 
with energy between E and E+dE at time t. As time goes 
on, a fraction of these systems reaches the energy Eb and 
undergoes gravitational collapse towards Eq. Therefore, 
we adopt the boundary condition 



PiEB,t)^0. 



(132) 



Our aim, now, is to estimate the current of diffusion past 
Eb and the typical lifetime of metastable states. 

6.2. The stationary solutions 

The stationary solutions of Eq. (|131|l are of the form 



(133) 



where J < is the current of diffusion in energy space. 
Using the boundary condition H132|) . the solution of Eq. 
(tn^ reads 



P{E) = Je-^^(^) / 



Eb ^ldF{x) 



dx. 



E D{x) 

The current of diffusion can therefore be expressed as 



(134) 



J 



(135) 



IEa D(E) 

To estimate the probability Pa, we shall approximate the 
curve F{E) close to A by a parabole. We thus make the 
expansion 



(128) F{E) = F{E^_ 



1 



2TCa 



{E - Ea? 



(136) 
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Fig. 15. Illustration of the barrier of free energy corre- 
sponding to Fig. El Point A corresponds to the gaseous 
metastable state and point B is the unstable solution cre- 
ating a barrier for reaching the condensed state C. 



where we have used F"{E) = -Tf3'{E) = 1/TC where 
C = dE/dT is the specific heat. Therefore, 



Pa = 



1 



/3 



-f OO 



-(E-Ea 



dE 



(137) 



On the other hand, the integral in Eq. (|135|l is dominated 
by the value of the integrand close to B. Making the same 
quadratic expansion as in Eq. (|136|l . we get 



Eb ^I3F{E) 



Ea 



D{E) 



dE: 



I f,l3F{EB 

d{Eb) y_ 



e^'^^-^^'l'dE 



oPF{Eb) 



2D{Eb) 13 



(138) 



where we recall that Cb < for the unstable solution. We 
thus obtain the expression of the current 



n^CA\CB\ 



-Fa) 



This expression involving the barrier of free energy AF, is 
similar to the one obtained by Kramers (1940) in his clas- 
sical study. In our case, the parameters have a thermody- 
namical interpretation while Kramers considers a dynam- 
ical system in a potential U{x). 

6.3. The escape time 

The preceding approach assumes that the population of 
systems that are introduced at A is continuously renewed 
so as to counterbalance the population of systems that are 
lost at B and maintain a stationary regime. We shall now 



relax this simplifying assumption and look for decaying 
solutions of Eq. 113111 of the form 



P{E,t) 



-At 



h{E), 



where h{E) satisfies the differential equation 

= -\h. 



_d_ 

dE 



D[^+PHE'iE) 



(140) 



(141) 



Assuming that close to Ea the system is at equilibrium, 
the foregoing equation can be integrated into 



dh , , , A 

dE+^'^^'^^-'lW) 



h{E')dE'. 



(142) 



Ea 



In usual situations, the eigenvalue A is expected to be 
small as it corresponds to the inverse lifetime of the 
metastable states. We thus consider the perturbative ex- 
pansion of h in powers of A and write 



h = ha -\- Xhi + 



(143) 



Substituting this expansion in Eq. 1142() and identifying 
terms of equal order, we obtain the differential equations 



+ PhoF'{E) = 0, 



dho 
llE 



Ea 



ho{E')dE'. 



The first equation integrates into 
Substituting this result in Eq. H145|) . we obtain 



dhi 

Ie 



-(3h,F'{E) = - 



K 



D{E) 



(144) 
(145) 

(146) 

(147) 



Ea 



The solution of this first order differential equation can be 
written 



(139) x'iE) 



where the function x is defined by 



D{E) 



(148) 



(149) 



Ea 



with x{Ea) — 0. Therefore, in the approximation A <C 1, 
the solution of Eq. H131() is 



P{E, t) = i^e-^*e-^^(^' [1 - XxiE)] . 



(150) 



The eigenvalue A is determined by the boundary condition 
(fn^ yielding 



A 



1 



(151) 



xiEsy 

Therefore, the lifetime of the metastable state is given by 

tufe^X{EB). (152) 
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We can now try to simplify this expression. First, we ap- 
proximate Eq. (|149|l by 



X'{E) 



D{E) y_ 



-0F{EA)g^(E-E_. 



dE 



(153) 



P D{E) • 

After integration, we obtain 

With the additional approximation 



^TTr^^~^^^'^"^'dE, (155) 



DiEe) 



we finally get 
A 



P^D{Eb) -P(Fb-Fa) 



(156) 



Therefore, the lifetime of a metastable state behaves as 



AJ 



(157) 



We note that the expression of A is similar to the expres- 
sion (|139|l obtained for the current J. The connexion is 
the following. In the non-stationary case, the current of 
diffusion at Eb is 

Jb^-D{Eb)(^^+PPF'{E)^ 



Eb 



= -\e-^'DiEB) { + l3h^F'{E) 



= Ae--^*/ ho{E')dE'c^X P{E',t)dE'. (158) 
Jea Jea 

Hence, normalizing the current by the exponential decay 
of the density probability, we get 



Jb 



(159) 



which is equivalent to Eq. p39|) . 

In the preceding analysis, we have worked in the canon- 
ical ensemble because the Brownian model (|124ll - (|125f) is 
easier to study than the A^-stars Hamiltonian model, while 
exhibiting qualitatively the same phenomena (phase tran- 
sitions, metastable states etc.). We expect to have sym- 
metric expressions in the microcanonical ensemble with 
the correspondance E ^ j3 and S ^ J . This study is left 
for a future work. 



7. Conclusion 

In this paper, we have completed previous investigations 
concerning the statistical mechanics of self-gravitating 
systems in microcanonical and canonical ensembles. The 
microcanonical ensemble is the proper description of 
isolated Hamiltonian systems such as globular clusters 
(Binney & Tremaine 1987). The canonical ensemble is 
relevant for systems in contact with a heat bath of 
non-gravitational origin. It is also the proper description 
of stochastically forced systems such as self-gravitating 
Brownian particles (Chavanis, Rosier & Sire 2002). We 
have justified the mean-field approximation, in a proper 
thermodynamic limit N — > +oo with rj — f3GMm/ R and 
e = ER/GM"^ fixed, from the equihbrium BBGKY hier- 
archy. In this thermodynamic limit, the equilibrium state 
is determined by a maximization problem: the maximiza- 
tion of entropy at fixed mass and energy in the micro- 
canonical ensemble and the minimization of free energy 
at fixed mass and temperature in the canonical ensemble. 
This determines the most probable macroscopic distribu- 
tion of particles at equilibrium. This can also be seen as a 
saddle point approximation in the functional integral for- 
mulation of the density of states and partition function. 
We have shown that the saddle point approximation is 
less and less accurate close to the transition point since 
the condition N\E ~ Et\ ^ 1 (in microcanonical ensem- 
ble) or A^|r — Ttl ^ 1 (in canonical ensemble) must be 
satisfied. 

We have also argued that the lifetime of metastable 
states (local entropy maxima) scales as exp(iV) due to the 
long-range nature of the interaction. Therefore, the im- 
portance of these metastable states is considerable and 
they cannot be simply ignored. Metastable states are in 
fact stable and they correspond to observed structures 
in the universe such as globular clusters. The preced- 
ing estimate must, however, be revised close to the criti- 
cal point. By solving a Fokker-Planck equation, we have 
shown the the lifetime of metastable states is given by 
the Kramers formula involving the barrier of entropy or 
free energy. These barriers have been calculated exactly 
close to the Antonov energy Ec (in microcanonical en- 
semble) and close to the Jeans-Emden temperature Tc 
(in canonical ensemble). We have obtained the estimates 
tiife ~ exp{1.726 A^(Ac — A)^/^} (in microcanonical en- 
semble) and tiife ~ exp{0.339 N{r]c — ?7)^^^} (in canonical 
ensemble) so that the lifetime decreases as we approach Ec 
or Tc- This implies that the collapse will take place slightly 
above E^ or at an energy A; = Ac(l - 2.077 N^"^/^) or 
temperature rji — ryc(l — 0.816 N~'^/'^). Similar conclusions 
have been reached by Katz & Okamoto (2000). Yet, these 
predictions do not seem to be consistent with the Monte 
Carlo simulations of de Vega & Sanchez (2002), although 
they find that the collapse indeed takes place slightly be- 
fore the critical point. Independent simulations are under 
preparation to check that point. 

Finally, a part of our discussion was devoted to an- 
swer the critics raised by Gross (2003,2004) in recent 
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comments. This author argues that the microcanonical 
entropy SmicroiE) and the microcanonical temperature 
PmicroiE) must be single valued. This is true in a strict 
sense, but the problem is richer than that because of the 
existence of long-lived metastable states. Therefore, the 
physical caloric curve/series of equilibria (3{E) is multi- 
valued and leads to "dinosaur's necks" and special "mi- 
crocanonical phase transitions" (Chavanis 2002). This is 
specific to systems with long-range interactions in view 
of the long lifetime of metastable states (local entropy 
maxima). These results have stimulated a general classi- 
fication of phase transitions by Bouchet & Barre (2004). 
Microcanonical phase transitions (as in Fig. ^ have not 
been fully appreciated by Gross and his collaborators be- 
cause their studies (e.g., Votyakov et al. 2002) consider a 
large small-scale cut-off for which the caloric curve looks 
like Fig. |21 and is univalued. If these authors reduce their 
small-scale cut-offs, they will see "dinosaurs" appear! 
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Appendix A: Justification of the mean-field 
approximation from the equilibrium BBGKY 
hierarchy 



As is well-known, each equation of the hierarchy involves the 
next order distribution function. We now decompose the two- 
body and three-body distribution functions in the suggestive 
forms 



P2 (r 1 , ra ) = Pi (r 1 )Pi (ra ) Pa (r 1 , ra ) , 



(A.6) 



Pa (ri , ra, rg) = Pi (ri)Pi (r2)Pi (rg) P^(ri , r2)Pi (rs) 
+P2 (r 1 , ra ) Pi (ra ) + Pa (ra , rg ) Pi (r i ) + P^ (r i , ra , rg ) , (A. 7) 

where are the cumulants. We shall consider the thermody- 
namic limit A'^ — > -l-oo with fixed 77 — fiGNm^ / R. In this limit, 
it can be shown that the non trivial correlations P^ are of or- 
der iV-("-i). Here, we shall just establish this result for the 
two-body distribution function Pa . Substituting the decompo- 
sitions ljA.6fl and l|A.7^ in Eqs. I|A.4^ and l|A.5fl and assuming 
that P3 is negligible (this corresponds to the Kirkwood ap- 
proximation in plasma physics) the equation for the two-body 
distribution function becomes after simplification 

— (ri,ra) = -/3Pi(ri)Pi(r2)^ 

-/3P^(ri,ra)^-/3iVP^(ri,ra) / Pi(r3)^d^r3 
OTi J an 

-/3AfPi(ri) y"p^(ra,r3)^d'r3. (A.8) 

We thus find that Pi ~ 1 and Pa ~ /?u ~ PGm'^/R = rj/N = 
0{1/N). Therefore, in the limit N — > -l-oo, the two-body dis- 
tribution function is the product of two one-body distribution 
functions: 



In this Appendix, we show that the mean-field approxima- 
tion is exact for self-gravitating systems in a properly defined 
thermodynamic limit N +00 with rj — f3GMm/ R and 
A = —ER/GM^ fixed. In the canonical ensemble, the equi- 
librium A''-body distribution function is given by 

p _ 1 -/3i7(ri,...,rjv) If. -i^, 

where we only consider the configurational part (the velocity 
part, which is just a product of Maxwellians, is trivial). Here, 
C/(ri, rjv) = X]i<j where Uij = -G/|ri -Vj\ is the grav- 
itational potential (we can also use a soften potential in order 
to regularize the partition function). Taking the derivative of 
Eq. lA.H with respect to ri, we get 



dPN _ du 
"a — ~ -p-nv — • 

OTi OTl 



(A.2) 



From this relation we can obtain the full equilibrium BBGKY 
hierarchy for the reduced distribution functions (Chavanis 
2004b). Restricting ourselves to the one and two-body distri- 
bution functions 



P,(ri, Fj) = / Pjv(ri, .... r]v)d^rj+i...d^rjv. 



with j = 1, 2, we find 

aPi 
9ri 



— i(ri) = -/3(iV-l) 



P3(ri,ra)-;^ — d ra. 



9P2 , ^ „ „ . ^ dui2 
-75— (ri,ra = -/3Pa(ri, ra)^; — 
ari dv\ 

-l3{N-2) j P3(ri,ra,r3)^d^3. 



(A.3) 



(A.4) 



(A.5) 



Pa(ri,ra) =Pi(ri)Pi(ra)-HO(l/Ar). 



(A.9) 



This justifies the exactness of the mean-field approximation 
for self-gravitating systems. Note that rj/N can be identified 
as the gravitational version of the "plasma parameter". This 
is similar to the remark of Lundgren & Pointin (1977) for the 
point vortex gas. Now, plugging this result in Eq. I|A.4^ . we 
find that, for A*' ^ -f-oo. 



|^(ri) = -/3APi(ri) /Pi(ra)^d^ra. 

OTl / OTl 



(A.IO) 



Integrating with respect to ri and introducing the mean den- 
sity p(r) = . mSijCi — r)) = NmP\{v), we obtain the 
Boltzmann distribution 



Ae 



(A.ll) 

where "I'(r) = J p{r')u{r — r')d^r' is the self-consistent gravita- 
tional potential. Adding the gaussian velocity factor, we obtain 
the Maxwell-Boltzmann distribution 



f = A'e 



(A.12) 



As we have seen in Sec. El the distribution function ljA.12|l 
can also be obtained by minimizing the Boltzmann free energy 
Pb[/] at fixed mass M and temperature T. This method pro- 
vides a condition of thermodynamical stability S^F > 0, which 
is not captured by the equilibrium BBGKY hierarchy. To get 
the condition of stability, we need to consider time-dependent 
solutions, i.e. the non-equilibrium BBGKY hierarchy. Indeed, 
the thermodynamical stability is related to the dynamical sta- 
bility with respect to the Fokker-Planck equation (Chavanis 
2004c). 
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The equation for the two-body distribution function lA.SII 
is complicated because the one-body distribution function is 
spatially inhomogeneous. It may be of interest, however, to 
advocate the Jeans swindle and consider, formally, the case of 
an infinite homogeneous self-gravitating system (this can be 
made rigorous in a cosmological context; see Kandrup 1983). 
Making the drastic approximation Pi = p/M where p is a 
constant, Eq. HA. St simplifies into 



dh 
dri 



(ri,r2) 



dui2 „p f dui3 s 

OTi m / OTi 



(A.13) 



where the second term in the right hand side of lA.SII . of order 
I/N'^ , has been neglected. The correlation function h is defined 

by 



(A.14) 



where x = ri — r2. Taking the divergence of Eq. ljA.13^ and 
using Alt = 47rG'm^(5(ri — r2), we obtain 



A/i + k% = --47rG/3m^5(x), 



(A.15) 



where kj = (47rGm/3p)^'''^ is the inverse of the Jeans length. 
This equation is easily integrated to yield 



h{x) = PGm 



2Cos{kjx) 



(A.16) 



This is the counterpart of the Debye-Hiickel result in the grav- 
itational case (Kandrup 1983). We emphasize that the above 
results are valid for other systems with long-range interactions 
(Chavanis 2004b). In particular, for the HMF model for which 
a homogeneous phase rigorously exists, we find by the same 
method that h{6) — -^j^^y^ cos6', where Pc = is the 
critical inverse temperature. In particular, the correlation func- 
tion diverges close to the critical point where the homogeneous 
phase becomes unstable, so that the mean-field approximation 
ceases to be valid. We expect a similar behavior for inhomoge- 
nous self-gravitating systems close to Tc. 

Considering now an isolated Hamiltonian system, the A^- 
body microcanonical distribution function is given by 



Pjv(ri,vi,...,rjv,vjv) 



giE) 



5{E - H{ri,\ri, rjv, vjv)). 



(A.17) 



From this expression it is easy to write the equilibrium BBGKY 
hierarchy (Chavanis 2004b). The first two equations of this 
hierarchy are 



9 Pi 



dui2 1 d 
dri giE)dE 



9{E)P2(l,2) 



d\2), 



(A.18) 



dP2 
dri 



(1,2) 



-(iV-2) 



9iii2 1 d 
' dri g{E) dE 

dui3 1 d 



g{E)P2(l,2) 
{E)P-i{l,2,3) 



dn g(E)dE 
where we have written (j) = (rj,Vj). Now, 



_1 d_ 

g{E) dE 



9{E)P, 



' ' ' + dE ' 



d^(3), (A.19) 



(A.20) 



where P = dS/dE and S{E) = lng{E). The ratio of dPj/dE 
on pPj is of order -g^ = -j^^- Therefore, in the thermody- 
namic limit A'' +00 with A, -q fixed, the second term in 
the r.h.s. of Eq. 



HA. 2011 is always negligible with respect to the 
first. To leading order in A'^, we obtain the same equations as in 
the canonical ensemble. Therefore, the mean-field approxima- 
tion is exact and leads to the Boltzmann distribution lA.llll . 
Observing that 



dPi 
9vi 



(1) 



1 



d 



g{E) dE 



[g(E)Pi{l)]vi, 



(A.21) 



and taking the 

g/3i,J/2 Combined with Eq 



oo limit, we find that Pi ~ 
lA.llll . this leads to the Maxwell- 
Boltzmann distribution IIA.12L Therefore, the equilibrium 
BBGKY hierarchy in the microcanonical ensemble leads to 
the same result 1A.12II as in the canonical ensemble. As indi- 
cated previously, the inequivalence of ensembles will appear by 
considering the non-equilibrium BBGKY hierarchy. The ther- 
modynamical stability in the microcanonical ensemble is con- 
nected to the dynamical stability with respect to the Landau 
equation (see Chavanis 2004c) which can be deduced from 
the non-equilibrium BBGKY hierarchy to order 1/A'^ (Balescu 
1963). 

Appendix B: Derivation of Eq. (|128|> 

The phase space hypervolume with energy less than E is de- 
fined by 



JV 

I{E) = I Yld^.d^Vi 

' H<E ,_ 1 



(B.l) 



Integrating over the velocities and using the fact that the ki- 
netic term in the Hamiltonian is quadratic, a standard calcu- 
lation yields 



(B.2) 



[E-U{ri,...,r^)\^\[d\,, 

i=i 

where A = (2/m)''^''^ Vsjv and Vn is the volume of a unit- 
hypersphere in a space of dimension n. The density of states 
g{E) — dl /dE is therefore given by 

g[E) = ^A I [E-U{ri,...,rr,)]'^~'l[d'r.. (B.3) 

i — l 

Assuming now that Pjv(ri , Vi, rjv, vjv, t) ~ PN{E,t), and 
substituting this ansatz in the A''-body Fokker-Planck equation 
Eq. 1126L we obtain after simplification 



+3Am(D^+^Pjvj, 



(B.4) 



where the term in bracket is X/i=i"?- note that Pjv = 
PN{E,i) is not an exact solution of 1126II . as expected. To get 
rid of the dependence in ri, tat, we shall average Eq. IIB.4t 
over the hypersurface of energy E using 

HE-U)'^-'X[ri,...,r.-E)Y{l^d-r.. 

\{E-U)^ n.^i'^r. 
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according to Eq. IIB.311 . This gives 
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+ZMg{E)[D^+^P^ 



(B.6) 



Using g{E) — dl / dE{E), we can put this equation is the form 

Appendix C: Rotating self-gravitating systems 
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Sec. 12.31 the conservation of angular momentum is equivalent 
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AL/'] = j fr'xVd\'d\'. (C.l) 

Now, repeating the argumentation of Sec. |2l the density of 
states 

g{E, L) = / 5{E - H)5 - ^ mr, x v.^j J| (fv,d\„ 
•' ^ j=i ^ i=i 



(C.2) 



can be written formally as 

5(£,L)= j Vf e'^f^S{E-E[f]) 
xS{M~M[f])5{L~L[f]). 
Similarly, the partition function 

/JV 
^-im+pn-Y,^^^ ™r. XV. -Q ^^^^^^ 
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can be written as 

Z{P,n)^ I e^'^''+''"-^g{E,-L)dEd^-L, 



(C.5) 



Z{/3, f2) = / Vf e-''^i(5(M - M[f]) 



(C.6) 



where J[f] = S[f]- l3E[f] + /3n-L[f] is the free energy. In order 
to apply the saddle point approximation, we just need to im- 
pose that A = -ER/GM'^, r] = f3GA4m/R, A = L/iGM^R^^^ 
and uj = fl(R'^ /GMy^^ remain of order unity in the limit 
A'' +00 (in the case of self-gravitating fermions, we also 
need to impose that fi = r/oVXPAIR^ is fixed and in the case 
of a soften potential that e = ro/R is fixed). This defines 
the thermodynamic limit for rotating self-gravitating systems. 
The corresponding scalings are given in Chavanis & Rieutord 
(2003). In particular, 5 ~ iV and J ^ N. Therefore, in the 
N ^ oo limit, we have to maximize S[f] at fixed E, M and 
L in the microcanonical ensemble and we have to maximize 
J[f] = S[f] - l3E[f] + p n- L[f] at fixed (3, M and fl in the 
canonical ensemble. Computation of rotating self-gravitating 
systems in relation with statistical mechanics have been per- 
formed by Votyakov et al. (2002) for a classical gas on a lattice 
and by Chavanis & Rieutord (2003) for fermions. 



